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ABSTRACT 

A new approach to the cosmological recombination problem is presented, which completes 
our previous analysis on the effects of two-photon processes during the epoch of cosmolog- 
ical hydrogen recombination, accounting for ns-ls and nd-ls Raman events and two-photon 
transitions from levels with n > 2. The recombination problem for hydrogen is described us- 
ing an effective 400-shell multi-level approach, to which we subsequently add all important 
recombination corrections discussed in the literature thus far. We explicitly solve the radia- 
tive transfer equation of the Lyman-series photon field to obtain the required modifications 
to the rate equations of the resolved levels. In agreement with earlier computations we find 
that 2s- Is Raman scattering leads to a delay in recombination by AN e /N e ~ 0.9% at z ~ 920. 
Two-photon decay and Raman scattering from higher levels (n > 3) result in small additional 
modifications, and precise results can be obtained when including their effect for the first 3-5 
shells. This work is a major step towards a new cosmological recombination code (CosmoRec) 
that supersedes the physical model included in Recfast, and which, owing to its short runtime, 
can be used in the analysis of future CMB data from the Planck Surveyor. 

Key words: Cosmic Microwave Background: cosmological recombination, temperature 
anisotropies, radiative transfer 



1 INTRODUCTION 

The Planck Surveyo£] is currently observing the temperature and 
polarization anisotropies of the cosmic microwave background 
(CMB), and scientists all over the world eagerly await its first data 
release, scheduled for early 2011. With Planck data sets cosmol- 
ogists will be able to determine key cosmological parameters with 
unprecedented precision, making it possi ble to distinguish betw een 
the various models of inflation (e.g. see IKomatsu et al J 1 20 id , for 
recent constraints). 

Ov er the past 5 yea rs, v ario us groups (e.g. 
see iDubrovich & Grachev! 120051: IChluba & Sunvaev 



Kholupenko & Ivanchild 
Wong & Scottl 120071: 



20061: 



Switzer & Hirata 
Rubino-Martfn et alj 120081: 



Karshenboim & Ivanovl 120081: iHiratal 120081: ICh luba & Sunvaev 
2008 ; Jentschural l200g|; lLabzowskv et al.l 120091 : iGrin & Hirate 
20 id : lAli-Haimoud & Hiratd l2oTrF have realized that the un- 
certainties in the theoretical treatment of the cosmological 
recombination process could have important consequences for 
the analysis of the CMB data from the Planck Surveyor. It was 
shown that in particular our ability to measure the precise value 
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of the spectral index of scalar perturbations, n s , and the baryon 
content of our Universe will be com promised if modification s to 
the recom bination model by Recfast dSeager et al J 19991 . 120001) are 
neglected jRubino-Martm et ai1l2010h . 

To ensure that uncertainties in the cosmological recombina- 
tion model do not undermine the science return of the Planck 
satellite, it is crucial to incorporate all important processes lead- 
ing to changes in the free electron fraction close to the maxima 
of the Thomson visibility function dSunvaev & ZeldovichI [l97bT) 
by more than ~ 0.1% into one recombination module. The main 
obstacle towards accomplishing this so far was that detailed re- 
combination calculations took too long to allow accounting for 
the full cosmological dependence of the recombination corrections 
on a model-by-model basis. Th is led to the introduction of im- 
prove d fudge factors to Recfast dWong & Scottll2007l ; Wong et al 



20081) . or multi-dimensional interpolation schemes dFendt et al 



20091 : 1 Rubifio-Martm et ai]|2010l) . that allow fast but approximate 



representation of the full recombination code. 

Although it was already argued that for the stringent error-bars 
of toda y's cosmological parameter s such approaches should be suf- 
ficient dRubifio-Martfn et~al]|2010h . from a physical stand point it 
would be much more satisfying to have a full representation of the 
recombination problem, that does not suffer from the limitations 
mentioned above, while capturing all the important physical pro- 
cesses simultaneously. Furthermore, such a recombination module 
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increases the flexibility, and allows us to provide extensions, e.g., 
to account for the effect of dark ma tter annihilation, energy injec- 



IPadmanabhan & Finkbeinei] 20051; Chlubd 


2010 


), or the variation 


of fundamental constants ( Kaplinghat et al. 


1995 


tlGalli et al.l2009l: 



In this paper, we describe our new approach to the recombina- 
tion problem, which enables us to fulfill this ambition by overcom- 
ing the problems mentioned above. Our code, called CosmoRecQ, 
runs in 1 - 2 minutes for a given set of cosmological parameters as 
it stands and can be optimized further to run well below a minute, 
eliminating the need for fudge factors to solve the recombination 
problem. One of the key ingredients that facilitates this increase 
in speed is the effective multi-level appro ach, which was proposed 
recently bv lAli-Haimoud & Hirat3 ( bold) . 

We also extend our previous analysis on the effects of two- 
photon pr ocesses during the cosmological re combination epoch of 



hydrogen dChluba & Sunvaevl 120081 201 Obi) to account for ras-ls 
and nd- Is Raman scattering and two-photon transitions from levels 
with n > 2. The radiative transfer equation for the Lyman-series 
photons during hydrogen recombination is solved in detail using 
a PDE solver that we developed for this purpose and can accom- 
modate non-uniform grids (see Appendix |B1 for more details). Our 
results for the effect of Raman scattering on the recombinat ion dy- 
namic s are in good agreement with earlier computations dHiratal 
120081) . Furthermore, we show that two-photon decays from levels 
with n > 4 - 5 can be neglected and Raman scattering is only im- 
portant for the first few shells. 

The main difficulty with two-photon and Raman processes 
during the recombination epoch is the presence of resonances in 
the interaction cross-sections related to normal '1 + 1' photon tran- 
sitions th at are already included into the multi-level recombina- 
tion code jffirata & Switzerl2008l ; lchiuba & Sunva^l2008l:lfflratal 



l2008l:IChluba & Sunvaevll2010bh . Unlike for the 2s- Is two-photon 
decay, all the higher ns-ls and nd-ls two-photon channels include 
'1 + 1' photon sequences via energetically lower Lyman-series res- 
onances, i.e., ns/d «-» kp «-» Is with k < n. Similarly, for ns-ls and 
nd-ls Raman-scattering events all higher Lyman-series resonances, 
i.e., ns/d <-» &p <-> Is with k > n, appear. Therefore, special care 
has to be taken to avoid double-counting of these resonances in the 
rate equations of the multi-level atom, as we explain in § 13.4.31 
§ l3A4l § 13331 and § |3~5H 

In §[2] we outline our principle approach to the recombination 
problem. The terms for the radiative transfer equation that allow to 
take all important recombination corrections into account are de- 
rived in §[3] We then solve the evolution equation for the high fre- 
quency photon field during the recombination epoch, and illustrate 
the different changes in §|4] In §|5] we discuss the different correc- 
tions to the ionization history, and we present our conclusions and 
outlook in §[6] 



2 PERTURBATIVE APPROACH TO SOLVING THE 
FULL RECOMBINATION PROBLEM 

2.1 General aspects of the standard recombination problem 

The cosmological recombination problem consists of determining 
an accurate estimate of the free electron fraction, X s = N s /Nji, as a 



function of redshift. Because of particle conservation, and the num- 
ber of electrons in excited states of H i and He i being negligible, 
one may writ^] 



N c « JV H [1 - Xfj + Nulfyc ~ Xul 



(1) 



where Nn denotes the number density of hydrogen nuclei, and 
,/kc = ^Vhc/^h is the fraction of helium nuclei. The populations 
of the different levels are given by Xf = Nf/N H , where 'a' indicates 
the atomic species. Furthermore, Xj = Xf for convenience. 

Equation (T) implies that the recombination problem reduces 
to finding solutions to X\ . For hydrogen, the standard rate equation 
describin g the evolution of the grou nd-state population has the form 
(see alsoi lSeager et al .1 19991 l2000l) 



dX» 
~dT 



(2a) 
(2b) 



k 

A*£„i s = A 2sU [Xg - Xf s e- hv * lkT r] 

A^„ ls = + »£„) [xg - 3 Xf s n*p ls ] . (2c) 

Here h^u is the mean photon occupation number over the Lyman-fc 
line profile, A, the atomic rate coefficients for spontaneous emis- 
sion, and nj^ ls is the occupation number of the CMB blackbody 
photons at the Lyman-A' transition frequency = Vkpu- 

The solution of Eq. ®, depends on the level populations 
of the 2s- and Arp-states. In addition, the photon distribution in 
the vicinity of every Lyman-resonance has to be known, to de- 
fine n*pi s . n^pis is often estimated by the Sobolev approximation, 
which, however, breaks down during recombination, leading to 
non-negligible corrections to the recombin ation dynamics (e.g. see 
IChluba & Sunvaevll2009bl . l2010bl l2009ij) . The rate equations for 
the 2s- and £p-states themselves can, in principle, be explicitly 
given. But here it is only important to realize that these lead to a 
large network of rate equations which depends on the populations 
of all other excited levels. To complicate matters further the elec- 
tron temperature, T c , enters the whole problem via recombination 
coefficients, ctj{T c , T y ), to each level where T y is the photon tem- 
perature. 

The evolution of T c is described by one simple differential 
equation, which accounts for the cooling of electrons caused by the 
Hubble-expansion, and the energy exchange with CMB photons via 
Compton scattering. Othe r processes, e.g., su ch as Bremsstrahlung 
cooling, are subdominant dSeager et alj 2000). 



2.1.1 The effective multi-level approach 

Recently, lAli-Hafmoud & Hiratal i20ld) suggested to simplify the 
recombination problem to a subset of levels that need to be fol- 
lowed explicitly. Here we shall call the members of this subset 're- 
solved' levels. This approach enables us to account for the effect 
of recombinations due to highly excited states (n > 100), without 
actually solving for all these level populations explicitly. The ratio- 
nale being that except for the optically thick Lyman-series transi- 
tions, all other rates are mediated by the CMB blackbody photons, 
and hence only depend on the photon and electron temperatures. 

The downside of this simplification entails the need to tab- 
ulate effective rate coefficients as a function of T y and T e prior 
to the computation. This however needs to be done only once, 



This code will be available at www.Chluba.de/CosmoRec. 



3 We assume to be at redshifts z <K 6000, well after doubly ionized helium 
recombines. 
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and given that the number of resolved states necessary for con- 
verging solutions is small, this effective multi-level approach re- 
sults in tremendous speed-up s for recombination calculations (see 
lAli-Haimoud & ffiratal lioTol . for more details). For this work we 
also implemented such an effective rate approach. The rate co- 
efficients for an effective 400-shell hydrogen ato m were com- 
puted using our most recent recombination code dChluba et al.l 
while helium is described with a multi-level treatment 
dChluba & Sunvaevll2010j) . 

Within this framework the choice of the number of resolved 
states depends on the extra physics that one intends to include. For 
example, in a minimal model for the hydrogen recombination prob- 
lem one should explicitly solve for th^| Is, 2s, and 2p level popu- 
lations in tandem with the electron temperature, T e . 

This minimal choice already allows us to include processes 
that affect the net rates in the 2s- Is two-photon channel and the 
2p-ls Lyman-or re sonance, e.g., the effect of stimulated 2s- Is 
two-photon decay dChluba & Sunvaev 2006b), o r the feedback o f 
Lyman-o- photons on the ls-2s rate l Kholupenko & Ivanchikl200^) . 
However, since we restricted ourselves to the Is, 2s and 2p states, 
corrections due to Lyman-/? or higher resonance feedback cannot 
be modelled as these wou ld require resolving np states with n > 2 
dChluba & Sunvaevl2007l) . We will return to these points in § 12.2.11 

2.1.2 Sobolev approximation for hk v { s 

In a multi-level approach the Sobolev approximation is invoked to 
obtain a solution for the photon-field around every resonance ap- 
pearing in Eq. d2"ct . The photon occupation number around each 
line is then given b)[f] 

An k v = Anl[l-e< 11 -^], (3) 

where An£ w «^ - nf pU , T k s is the Sobolev optical depth in the 
Lyman-fc resonance, and^' = (py(v)dv. Here (fiy(v) is the Voigt 
profile corresponding to a resonance, and the line occupation num- 
ber, w£, is defined as: 



2.2 Beyond the standard rate equation for Is 



3iV is 



(4) 



Consequently, a simple approximation for the mean occupation 
number is 

n*pi s = nl - P%Anl, (5) 

with P k & = [1 - e~ rs ]/Tg being the Sobolev escape probability. 

For the Lyman-a resonance Eq. ([3} results in a photon 

distribution that is rather unphysical (e.g. see discussion in 
I II k 

IChluba & S unvaev 2009b). This is primarily due to the assumption 
that every interaction with the resonance leads to a complete redis- 
tribution of photons over the whole line profile, which for typical 
values of t s during recombination couples the photon distribution 
from the line center up to frequencies in the Lyman-continuum. For 
conditions present in our Universe, photon redistribution over fre- 
quency is much less effective, most notably in the distant wings. 
Thus, it is important to distinguish between scattering, real emis- 
sion and absorption events, as we will discuss in §[3] 



4 Because the 2s and 2p states are usually close to full statistical equilib- 
rium, one could also eliminate either of these states using Xi p = 3X2 S and, 
as a result, closely resemble the normal Recfast code, now without requir- 
ing a fudge factor. 

5 We assumed that as v — » oo the distortion An„ — » 0, and that (as usual) 
the factor (vh/v) 2 ~ 1. 



As mentioned in § 12.1.11 within the effective multi-level approach 
the choice for the resolved states depends on the physics to be mod- 
elled in detail. For example, in order to include the full effect of 
Lyman-series feedback, say up to n = 8, the minimal model that 
follows Is, 2s, and 2p would at least have to be extended by all 
fcp-states up to 8p. 

Also, the inclusion of two-photon processes from higher levels 
and Raman-scattering, requires us to re-write equation l(2j in a more 
generalized form as, 



~dT 



Is _ Xj s .Rls- 



(6) 



where R^j are the rates between the levels i and j. These rates 
depend on atomic physics, the CMB blackbody, the electron tem- 
perature, and the solution for the Lyman-series spectral distortion 
introduced by the recombination process. 

To include two-photon corrections to the Lyman-series up to 
n < /2 milx , the important levels to follow are all the nd and ns-states 
with 2 < n < n max - The corresponding partial rates to the rcp-states 
drop out of the equations, and the Lyman-series emission and ab- 
sorption profiles, usually given by a Voigt function, will be replaced 
by the two-photon profiles for the ns «-> Is and nd «-> Is processes, 
and similarly for the Raman process. We will specify these correc- 
tions more precisely in the following sections. 



2.2.1 Accounting for corrections from radiative transfer effects 

Changes in the level populations, electron temperature and free 
electron fraction remain small (~ 1%), when different physical 
processes, which were ne glected in earlier treatments (e.g. see 
iRubino-Martm et al]|2010l for overview) are included. This justi- 
fies treating corrections to T c and the populations of resolved lev- 
els, Xj, as small perturbations. On the other hand, it has been shown 
that the changes in the pho ton field caused b y time-dependence 
dChluba & Sunvaev! l2009bh. line scattering dChluba & S unvaev 
| 2009a | ; Iffirata & Forbesll2009l), or two-photon corrections iHiratal 
l2008l : lchluba & Sunvaevll2010bl) . are non-perturbative. 

In i}3] we derive in detail the different correction terms for the 
photon diffusion equation and provide modifications to the net rates 
of the effective multi-level atom. The idea is to first solve the re- 
combination history using the effective multi-level approach in the 
'1 + 1' photon description, i.e. equate R^j = Rl*\ in Eq. {3), and 
then compute the solution to the photon field using the radiative 
transfer equation. This then leads to corrections in the net rates, 
which are used in computing changes to the recombination dynam- 
ics, and hence modify Eq. l[6). These corrections being small, de- 
mand only one iteration to converge. Detailed descriptions to the 
notations in the foll owing sections and part of the methods u sed 
can be also found in lChluba & Sunvaev! d2009bl. 1201 Obi, l2009al) . 



3 EQUATION FOR THE PHOTON FIELD EVOLUTION 
AND CORRECTIONS TO THE EFFECTIVE 
MULTI-LEVEL ATOM 

To account for all corrections to the cosmological recombination 
problem, it is important to follow the evolution of non-thermal pho- 
tons in the Lyman-series, which are produced during recombina- 
tion. These photons interact strongly with neutral hydrogen atoms 
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throughout the entire epoch of recombination, and their rate of es- 
cape from the Lyman-resonances is one of the key ingredients in 
accurately solving the recombination problem. 

The partial differential eq uation governing the evolut ion of the 
photon field has the form (see lChluba & Sun vaev 2009jl f° r a de- 
tailed discussion) 



1 


dAn v 


dAn,, 








-Hv—1 




c 


dt 


v OV 


/ 



C[An v ] em/abs + C[An v 



(7) 



where An v = j^AN v is the distortion in the photon occupation 
number, and a distinction is made between the collision terms 
leading to emission and absorption, C[A/3,,] cm / a bs, and scattering, 
C[A« v ] scat . As an example, the first term on the right hand side of 
the equation can account for two-photon corrections to the line pro- 
files, while the second term, electron and/or resonance scattering. 
The second term on the left hand side describes the redshifting of 
photons due to Hubble expansion, and plays a crucial role in the es- 
cape of photons from the optically thick Lyman-series resonances. 

In Eq. 10 the CMB blackbody has been subtracted, i.e., An v = 
n v - ;z pl , where rc pl = [e l,y/kT ^ - l] -1 is the blackbody occupation 
number, because the left hand side directly vanishes for a black- 
body with temperature T 7 (z) «. (1 + z). Also, spectral distortions 
created by Compton scattering off electrons with T c + T y will be 
extremely small for conditions in our Universal at least if there is 
no additional energy release. 

By changing the time-variable to redshift z using dz/ dt = 
-H{\+z), and scaling to dimensionless frequency x = v/vji, Eq. ^} 
reads 



dAn x 



dz 



dAn, 



[1 + z] dx 



- A ; {C[An 



vjem/abs 



+ C[A«„] S , 



(8) 



We will now discuss the 



where An x = V21 An,, and A ; 

terms describing the resonance and electron scattering. In § |3.1| we 
specify the different emission and absorption terms, which then in 
§|4]and ij5]are used to compute the corrections to the Lyman-series 
distortion and ionization history. 



where x is the dimensionless frequency and = jfr- « 40 . 
The first term on the right hand side describes photon diffusion and 
the second accounts for the recoil effect. 



3.1.1 Electron scattering 

Th e diffusion coefficient in the case of electron scattering i 
sec IZeldovich & Sunvaevll 19691 : Is azonov & Sunvaevll2000l) 



(TjNzC 

H[l+z] 



(10) 



where cr T » 6.65 x 10 25 cm 2 is the Thomson cross section. 
IChluba & Sunvaev! d2009al) pointed out that electron scattering has 
an effect only at the early stages of recombination (z > 1400). How- 
ever it is easy to include, and also has the advantage of stabiliz- 
ing the numerical treatment by damping small scale fluctuations of 
the photon occupation number caused by numerical errors, even in 
places where line scattering is already negligible. 

As can be seen from the form of the diffusion coefficient in 
Eq. i ll 01 ), the efficiency of electron scattering to a large extent is 
achromatic. This is in stark contrast to the case of resonance scat- 
tering, which is most efficient only in a very narrow range around 
the line center (see next paragraph). Furthermore, the number of 
free electrons drops rapidly towards the end of recombination, such 
that the Fokker-Planck app roximation is expected to break down 
dChluba & Sunvaevl2009al) . Nevertheless, the diffusion approxima- 
tion remajnssufficiej^OTCcaiiputations of the free electron fraction 
(see lAli-Haimoud et alfeOld) . 



3.1.2 Resonance scattering 

For resonance s cattering by a Lyman-fc line the diffusion coeffi- 
cient i s (e.g. see lBaskd[l^78bMlRvbicldl2006l : lchluba & Sunvaevl 
l2009al and reference therein) 



23*00 « P 



ofNuc 
H[l+z] 



kT e 



ma<r 



v 2 



-x <pf(x), 



(11) 



3.1 Inclusion of partial redistribution by resonance and 
electron scattering 

Here we provide the terms for the Boltzmann equation de- 
scribing the effect of (partial) photon redistribution by reso- 
nance and electron scattering. The form of the collision term 
for these cases wi thin a Fokker-Planck formulation was di s- 
cussed earlier (e.g. IZeldovich & Sunvaevl [19691: iBaskol [l97 8b a ; 
Rvbicki & dell' Antonioll 19941: 1 S azonov & Sunvaevll200Cl: iRvbickil 
20061 : IChluba & Sunvaevll2009d)~ Since we are only following the 
evolution of the distortion from a blackbody, and since it is clear 
that induced effects are neg ligiblfl one can readily write 



A ; C[An,,]U tt «i|-D(x) 



An x + ^(z)An x 



(9) 



The cooling of CMB photons by losing energy to keep electrons at T e ~ 
T y should lead to a y-distortion with y-parameter y ~ 10~ 10 - 10 -9 . The 
dissipation of energy by acoustic waves should lead to y ~ 10~ 8 . These can 
be neglected for our purpose. 

7 Eliminating the dominant term of the CMB blackbody leaves us with a 
term that is tiny because we are always in the distant Wien tail of the CMB 
at all times during recombination. 



where crf' p = 3 JIl ^ a,*p 



^ and Av D denote the resonant-scattering cross 
section and the Doppler width of the Lyman-fc resonance, respec- 



tively. For the Lyman-cc line <x, ~ 1.91 x 10 cm and Avj 



2P 



2.35 x 10~ 5 v 2] at z ~ 
is normalized as j o 

kp 

«5 



1 100. The Voigt profile / v p (x) = ip k V(x) Av* p , 



J^V* p (v)dv = 1. Where 
jCp = [v - v tp ]/Avp P is the distance to the line center in units of 
the Doppler width. 

The scattering probability of the Lyman-£ resonance, pj p , is 
determined by a weighted count of all possible ways out of the kp- 
state, Rk p (Ty), excluding the Lyman-series resonance being consid- 
ered, and then writing the branching ratio a^| 



Ah, 



A kp[s +R 



(12) 



yielding the probability for re-injection into the Lyman-fc reso- 
nance. The rates R^ p (T y ) and the scattering probabilities, p k S c(T y ), 
can be pre-computed, independent of the solutions obtained from 
the multi-level code. We detail the procedure b elow. 

Following iRvbicki & dell' Antonio! d 19941) . the diffusion coef- 
ficient is D oc 0y P (v). We neglect corrections due to non-resonant 



Stimulated emission for the Lyman-series has been neglected. 
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contributions (e.g. see lLeel2005l) in calculating the scattering cross 
section, which would lead to a different frequency dependence far 
away fro m the resonanc e (e.g. Rayleigh scattering in the distant 
red wing, |jacksonlll998l) . However, because it turns out that reso- 
nance scattering is only important in the vicinity of the Lyman-a 
resonance, this approximation suffices. 

It is also worth mentioning that Eq. dl II ) together with Eq. 10, 
in the limit of large optical deptlfl provides a viable description 
for the redistribution of photon over frequency. Unlike the case of 
complete redistribution (in which the reemission after each scatter- 
ing event occurs over the whole Voigt profile), in the recombina- 
tion epoch photons are only partially redistributed as a result of the 
Doppler motions o f the hydrogen at oms, so-called type-II redistri- 
bution as defined in iHurnmen fl%3) . 



the initial run of the effective multi-level recombination code. How- 
ever, the simplest way to define the ratio R^/R^ is to use the quasi- 
stationary approximation for the /zp-population (see details below). 
We note that in full thermodynamic equilibrium An*™ = 0, so that 
no distortion is created (An* = 0). 

Physically, Eq. l !14t includes two important aspects, which are 
not considered in the standard recombination calculation. Firstly, 
it allows for a distinction between scattering events on one side, 
and real emission and absorption events on the other. Secondly, 
it ensures conservation of blackbody spectrum in full thermody- 
namic equilibrium, even in the very distant wings of the lines. Refer 
IChluba & Sunvaevl fcOlObT) for a detailed explanation of the latter 
point, and on how this leads to one of the largest corrections in the 
case of Lyman-o- transport. 



3.1.3 Equilibrium solution for the scattering term 

Independent of the type of scattering, the equilibrium distribution 
with respect to the scattering term Eq. ((9} is given by 



An? 



An x (z)e 



-ffe)[A-A- ] 



(13) 



This is the expected Wien spectrum with the temperature defined 
by the electrons. The normalization An xo (z) is determined by the 
emission and absorption process. 

The optical depth to line scattering being extremely large in- 
side the Doppler cores of the Ly man-resonances (t s ~ 10 6 - 10 s 
during H i recombination) causes the photon distribution within the 
Doppler core to remain extremely close to equilibrium, An s ?' eq . 



3.2 Normal Lyman-fc emission and absorption terms 

In the normal '1 + 1' photon picture, the emission profile for each 
Lyman-series resonance is given by a Voigt-profile, tp k ^ , with Voigt- 
parameter a kp . Given the rate, Rt (T y , T c ), at which fres/f*°l electrons 
reach the fcp-state, and the probability of photon injection into the 
Lyman-& resonance, p^ m = p* p , the Lyman-fc line-em ission and ab- 
sorption term are (e.g. see lChluba & Sunvae v 2010b); 



dAn x 



dz 



Ly-* 



em/abs 



^ ~h[T7z] pl 



ffAnS. (14) 



The factor 1 /x 2 accounts for the translation from photon number 
to the occupation number because AN V oc ^An^, for which the 
Voigt-profile is defined. Also, p k ^ = 1 - p\f m is the death or the real 
absorption probability in the k' h Lyman-series resonance, and AnJjJ, 
and f kp are given by, 



A;;; 



Ap 



f? 



gkp R kp N ^ 

exp (li[v 



-hv a lkT y 



]/«V), 



(15a) 



(15b) 



where gnlgkp is the ratio of the statistical weights of the initial 



and final states. The function An^Ty, T c ) can in principle be pre- 
computed using the solution for the populations of the levels from 



9 During hydrogen recombination photons scatter efficiently off the 
Lyman-o; resonance out t o ~ 10 4 - 10 5 Doppler width (see Fig. 3 in 
IChluba & Sunvaevl l2009bh . However, the redi stribution of photons in th e 
distant damping wings still remains rather slow JChluba & Sunvae\l20093) . 

10 Electrons that did not enter the p-state via the Lyman-i: resonance. 



3.2.1 Computing An 1 ^ 

To solve the evolution of the photon field, one has to know at which 
rate photons are produced by the Lyman-resonance. This rate de- 
pends on A«en, as defined in Eq. dl5t . 

The rate equation for the evoluti on of the population in th e 
fcp-level has the form (see Appendix B lChluba~& Sunvaev 2010bh . 



dX kp 
dt 



dXkp 
dt 

Ly-k 



&Xkp 
dt 



Ly-k 



R kp X kp 



r kp 



I 



Ap 



— Atpi s Xi s J* p - A^pis Xkp l\? 

Skp 

f V k w P (v)e" lv -^ ]lkT rn v dv 
j ^(vXl+nJdv* I+jAvh)*!. 



(16a) 
(16b) 

(16c) 
(16d) 



In this picture the emission, absorption and resonance scattering 
terms are all treated simultaneously. In addition, the asymmetry 
between the emission and absorption profile in the Lyman-fc res- 
onance, as required by detailed balance, has been incorporated. 

Under quasi-stationarity, and using the definition of the death 
probability, p k J ', Eq. (16} yields 



gn R t P 

gkp Rkp x h 

such that with Eq. dl5t 



1 

~kp 

Pi 



Skp %kp _ j-kp 

gn X u 1 



+ J 1 



A;c 



Ap 



Ap 
Pi 



An? 



gkp Xkp 

gn~x7s~' 

kp 

1 + ?^P\ 

kp ^ 

p7 



r Ap 



+ I, 



-hvkp\slkT y 



(17) 



(18a) 



(18b) 



In the second step we used the normal Sobolev approximation, for 



which T 



Ap 



- P k 



' "tpi.J (f° r th e case °f Lyman-a compare 
also with Eq. (41) in fchluba & Sunvaevll2009bT) . 

From Eq. d 1 8bt we have Anem ~ An L , since for all Lyman- 
series resonances the second term in brackets is very small. Never- 
theless, for the total normalization of the line i ntensity close to the 
line ce nter, this small correction is important dChluba & Sunvaevl 
l2009bl) . in particular for the Lyman-a resonance. 

Also we would like to mention that for the Voigt parameter of 
the Lyman-A' profiles, a kp = A^ t /[4^Avp ], the total width of the line 
is used, where transitions induced by the CMB blackbody (e.g. to 
higher levels) are included. Numerically, it is possible to compute 
the total width for the Lyman-fc resonance with A k ^ t = A^^/p 1 ^- 
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3.3 The 2s-ls two-photon channel 

The 2s- Is two-photon channel provides the pathway for about 
60% of all electrons in hy drogen to settle into the ground state 
dChluba & Sunva ev 2006a). It therefore provides the most impor- 
tant channel in the cosmological recombination process. Here we 
treat the case of 2s- Is separately to illustrate the important approx- 
imations in the two-photon picture. The derivation outlined in this 
section is then used to obtain the corresponding terms for the two- 
photon processes from excited states with n > 2 (see § !3.4t . 

The net change of the number density of electrons in the 2s 
level via the 2s- Is two-photon channel is given by 



dX 2i -■ 



s X ls J <p 2 2 y s ii(v) n(v 2i - v) dv 
- A ti. x * J V%1 1 + [ 1 + "( y 2 1 - v)] dv, (19) 

where A^ ls = 8.2206 s -1 dLabzowskv etai]|2005h is the vacuum 
2s- Is two-photon decay rate, and tp"? denotes the 2s- Is two-photon 
decay profile normalized as J <fi^ dv = 1. Including all possible 
ways in and out of the 2s-level the net change of the number density 
of electrons in the 2s-state can be written as 



dX 2s 
dt 



dX 2s 
dt 



2)' 



A'-.,Y.-. 



(20) 



Here and R^ include the effect of all transitions to bound states 
with n > 2 and the continuum. 

In order to simplify the notation we now introduce 



(f(y)t= PV/(v 

Jo 

g\ ={nn' f 



)dv, 



^ = <[l+n][l+n'])f 



(21a) 

(21b) 
(21c) 



where /(v) is some arbitrary function of frequency and n = n(v) 
and ri = n(v') with v' = vn s - v. 

Then, under quasi-stationarity the solution for the population 
of the 2s-state is given by 



R 2s + ^2sls @2 



(22) 



In the multi-level approach the effect of stimulated two-photon 
emission is neglected leading to gf as 1. Also any CMB spec- 
tral distortion that is introduced by the recombination process 
(e.g. because of Lyman-o- emission) is omitted, implying g 2s « 

( n pi riP 1 ' \ * ~ exp(-hv2i/kT y ). In this approximation, the result 
from Eq. J22b becomes identical to the one obtained using Eq. l !2bt 
and Eq. l !20t . in the standard multi-level approach. 

However, in the recombination problem corrections to both 
gf and g 2i are important. For the stimulated two-photon emission 
only the occupation number given by the undistorted CMB black- 
body has to be considered and thus, 



2s,pl 



(23) 



which can be precomputed as a funct ion of temperature. Typical ly, 
g^" pl exceeds unity by a few percent dChluba & Sunvaevll2006bh . 

For Qrf one can make use of the fact that the distortions at 
either v or V are very small, so that 



n pl n pl ' + n pV An + n? 1 An' 



(24) 



Hence Eq. dl9l > can be re-written as, 



dX 2s 
~~dT 



A 2^ [x^" vli ' kTr - X 2s ] + A^ ls X u Ag? (25a) 
Agf = J ipll^'An + n^ An']dv =2 J''^ n pl ' Andv, (25b) 
where we defined the stimulated 2s- Is two-photon decay rate 

^L<?2 S ' Pl (Cf. 



within the CMB ambient radiation field as A^j* = _ _ 
IChluba & Sunvae v 2006b). Also Eq. l |25bt reflects the symmetry 
of the two-photon profile around v = v 2 i/2. 

Note that for g^' p only the CMB blackbody spectrum is im- 
portant and therefore can, in principle, be precomputed as a func- 
tion of photon temperature, T y . This also emphasizes the difference 
in the origin of the two terms of Eq. d25ab . g^' p> being the ther- 
mal contribution, while Agf arises solely because of non-thermal 
photons created in the recombination process. 

By comparing Eq. l |25t with Eq. j2b\ one can write down the 
correction to the 2s- Is net two-photon rate 

AKS^i. = A 2Ls A ^2 S ' Pl [Sue-**'** - X 2 s] + Aj Is X u Agf. (26) 



Here we introduced Ag^' p = gf' p '-l, which during recombination 
is of order ~ 1%. In Equation i25\ the integral Ag 2s depends on the 
spectral distortion introduced by the recombination process in the 
Wien's tail of the CMB blackbody. Including only the Lyman-a 
distortion provi des a manner in which to take its feedback effect 
into account (cf. iKholupenko & Ivanchikj|2006l) . 



3.3.1 The 2s-ls two-photon emission and absorption term 

In contrast to the Lyman-series channels, the terms for the photon 
radiative transfer equation in the case of the 2s- Is channel can be 
directly obtained from the net rate between the 2s and Is state as in 
Eq. d 1 9b , resulting in 



I ^1 

c dt 



^L^2r[l+»(v)][l+n(v 21 -v)] 

- A ^ln N i^ n ( v ) n ( y 2i - v). 



2)' 



(27) 



^2s = ' wnere the factor of two results from 
two photons being added to the photon field, and the 4n converts 
the units to per steradian. 

The reason for this simple connection to the net rate equation 
is related to the fact that every transition from the Is state to the 
2s level is expected to lead to a complete redistribution over the 
2s-ls two-photon profile. The main reason behind this assumption 
of redistribution is that the probability of coherent ls-2s scattering 
event is tiny because the 2s- Is decay rate is extremely small com- 
pared to the time it takes to excite a 2s-electron to higher levels or 
the continuum. 

However, some additional simplifications are possible. First, 
we can again replace the factors, [1 + n][l + »'], accounting for 
stimulated two-photon emission with those from the undistorted 
CMB blackbody. Furthermore, from Eq. l !24b . 



l/WV 



[1 +n pl ][l +ri 



pi'i 



An An' 

1 + — r + — 77 

npi , 7 pi' 



(28) 



Also, since the spectral distortions at very low frequencies are never 
important, one of the two terms in Eq. {28} (say the one related to 
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An') can always be omitted. Therefore we can rewrite Eq. J27t as 



1 dN v 



2y 



: 4L^2r[ A »em-/v S An v ] (29) 



where = + rc pl ][l + n p] '] and 

Ant = ^-e-'-->^ (30a) 

e -hv 2i IkTy 

f~ = 57 ; K ex P - v 2i]/ fcr r) • OOb) 
n p '(v) v ' 

If the term ^-jr is non-negligible, as might be the case at very low 
redshifts (z < 400), where the Lyman-a photons emitted at z ~ 
1400 redshifts into the 2s- Is absorption channel, one in addition 
has to subtract the term / v 2 , s An v / within the brackets of Eq. l |29t . In 
terms of x = v/v 2 \, z and An x the photon occupation number now 
evolves as, 



dz 



2)' 



(31) 



Here the 2s- Is cross section is given by o~^ ls = ' g 1 ^' 8 , and 2 J'* = 
4wv 2 i^'*. 

Eq. d31| > bears a striking resemblance to the emission and ab- 
sorption in the Lyman-series channels as in Eq. d!4t because one of 
the two photons that are involved in the 2s- Is two-photon process 
is drawn from the undistorted CMB blackbody spectrum, so that 
the evolution equation essentially becomes a one-photon equation. 
The difference is the absence of a death probability since practi- 
cally every electron that is excited to the 2s state will take a detour 
to higher levels or the continuum as p 2 ^ ~ 1 . 



3.4 Two-photon emission and absorption terms from excited 
levels with n > 2 

One of the most interesting modifications to the solution for the 
photon field is related to the deviations of the profiles for the 
different two-phot on emission and absorpti on channels from the 
Lorentzian shape dChluba & Sunvaevll2008h . For the recombina- 
tion problem only those one-photon sequences involving a Lyman- 
series resonance (e.g. 4d «-» 2p «-» Is) are important*"*! In this sec- 
tion we shall replace the standard '1 + 1' -photon terms for these 
channels with the full two-photon description that takes into ac- 
count the coherent nature of the procesQ 

W e generalize the approach detailed in IChluba & Sunvaevl 
feOlObl) for emission of photons close to the Lyman-Q' line to in- 
clude corrections around the Lyman-/? and higher resonances. 



11 All the other two-photon emission and absorption channel (e.g. 4d «-» 
2p «-» 2s) can be treated within a blackbody ambient radiation field, so 
that their net rate can be directly computed. Without deviations from the 
blackbody shape these will be extremely close to the normal '1 + 1' photon 
rates. Also they can only affect the net recombination rate as a 'correction 
to correction', because they only act on the electron 'feeding rates' into 
the main channels towards the ground state. A similar argument holds for 
Raman scattering events that do not directly connect to the ground state. 

12 Conditions persistent in the Universe at the recombination epoch makes 
collision s negligible, maintaining the coherenc e of the two-photon decay 
(e.g. see lChluba & Sunva ev 2008; Hirata 2 0o3) 



3.4.1 Net rates for two-photon transitions from excited s- and 
d-states 

The net change of the number density of electrons in the level 
j e {ns, nd} via the j- Is two-photon channel is given by 



dXj 
~~d7 



= ^-A%X u f <pf(v)n(v)n(v jl -v) 
-■VXx.jc.ir 



dv 



)[l+n(v)][l+«(v /1 -v)]dv. (32) 



Here ip'J denotes the profile for the j- Is two-photon decay, which 
can be computed as explained in Appendix lAl and is normalizecf"! 
as J ip 2y dv = 1. The (vacuum) two-photon decay rate is given by 



(33) 



The ratio of the statistical weights is gjlgu = 1 for the ns-states, 
and gj/gu = 5 for /id-states. Equation d33 1 > simply reflects the one- 
photon decay rates and branching ratios of all the '1 + 1' photon 
routes j — » np — > Is via intermediate p-states with n < nj. Stimu- 
lated emission induced by the CMB photons is not included in the 
definition of A~[ s , since it is taken into account differentially by the 
integrals in Eq. d32b . 

With notations defined in Eq. J2U . and following the proce- 
dure to derive Eq. j25\ , we can re-write Eq. d32t as 



dXj 



Si 



Xue -'"i iJkT y - x, 



mi 



J VjlJ2 



A/; dv, 



+ £L A 2 j y Xls Ag J l (34a) 
gu J 



(34b) 



where we defined the stimulated j-ls two-photon decay rate within 
the CMB ambient radiation field as A 1 ^* = Aj ls Q'f. 

The Q'f 1 depends crucially only on the CMB blackbody spec- 
trum and thus can be precomputed as a function of photon temper- 
ature, T r . On the other hand, like for the 2s- Is two-photon process 
(see Eq. {25)), AQ\ arises due to non-thermal photons, and hence 
depends directly on the solution for the photon field. 

In the normal '1 + 1' photon picture, the two-photon profiles 
can be considered as a sum of <5-functions and therefore 



.2*1+1),. 

A jls " Zj j k T> P ' 



(35) 



Here A* kp = Ajk p [l + n pl (v,x)], and the stimulated effect close to the 

Apis 



Lyman-series resonances has been neglected, i.e. 1 +n pl (vi- p i s ) ~ L 



3.4.2 Two-photon emission and absorption for excited s- and 
d-states 

The two-photon emission and absorption terms are obtained fol- 
lowing the steps in the derivation of Eq. J3 lb . For the j- Is two- 
photon channel one therefore obtains 



<9An, jU 

2y 



dz 



o- 2 :.N u c v 2 <jf y -' , 

■WTz\iM 



(36) 



13 Small correction to the normalization due to the two-photon description 
are neglected. 
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The j- Is two-photon cross section is given by <r.[ s = j^- -j^rf, 

and <^ y '* = An V/iy5 2r (v) [1 + M pl (v)][l + n pl (v')], where, because of 
energy conservation, V = vju - v. Also, 



gu 

e -hvj\lkT y 



(37a) 



„A[v-v ;1 ]/*J,, 



^ = ^ s TT^ wexp ( Mv -^)- (37b) 

Again we emphasize the resemblance of the equation above to that 
of the one-photon equation for the Lyman-series emission and ab- 
sorption channels as in Eq. d 1 4b - 



3.4.3 Correcting the Lyman-series emission and absorption 
terms in the radiative transfer equation 

Two-photon decays from a given initial state j e j«s, «dj involve 
Lyman-series resonances with k < n. For example, a 4d-ls two- 
photon emission event includes the effect of the Lyman-a and f} 
resonance. In the Lyman-series emission and absorption terms as 
in Eq. d 14t . these are already accounted for as ' 1 + 1' photon terms, 
when the profile is given by the normal Voigt function. 

To avoid the double counting of these transitions in the radia- 
tive transfer equation, two modifications are necessary: (i) all death 
probabilities, p kp , have to be reduced to account only for those 
channels that are not included in the two-photon description, and 
(ii) the Lyman-series emission rates have to be r educed for the same 
reason. This approach was also explained in IChluba & Sunvaevl 
d2010bh for the 3s- Is and 3d- Is two-photon process. Including only 
the j'-ls two-photon process (say for 3d- Is), the modified death 



probability and An^J, of the Lyman-A; resonance becomes, 



-kp 
Pi 



kp j,kp 



1 R + — 

L\ft e , 



3X ls R- - a> 

kp kp 
./'■'- P 



e-hvtylkTy 



(38a) 
(38b) 



where the partial death probability, p i , is given by 



of all '1 + 1' photon processes, j <-» np <-» Is, via Lyman-series 
resonances with n < iij. Double-counting can again be avoided by 
subtracting the corresponding '1 + 1' photon terms from the full 
j- Is two-photon rate. The remaining corrections can then be added 
to the effective multi-level code as additional rates which directly 
connects level j to the ground statj* 4 ! 

In the standard multi-level description of all _/'«-> np <-» Is 
sequences (iij > n), the contributi ons to the two-photon net ra te as 
in Eq. l |32| l, take the form (see also lChluba & Sunvaevll20fobT) 



~~AT 



dXj 
~At~ 



2r(l+l) 



I s,kp 



2*1+1) 



Skp y a J*f » 
Ai s Akpis p A ntpi 

gu 



6J v A 
A Is Aj kp p, 

gu 



x J A ~Jkpp 

kp » pl (V/*)«*pi 



kp 



Y-A' 

A J n jkp Psm 



v Ei 

At 



2y(l+U 



sip 



(40a) 



(40b) 



(40c) 



Equation J40al > is interpreted as electrons exiting level j via the 
route j — > np — > Is at a rate A*^ times the probability, pem (sec- 
ond term). Similarly, electrons reach state j from the ground state 
via the route Is — > np — > j, with the Lyman-A; excitation rate, 
Atpi s n^pis times the probability, , to then make the transition 
kp — » j (first term). Using Eq. d38c| ) and l |38e| l leads to Eq. l |40bt . 

Equation d40bt helps make the connection of the full two- 
photon net rate and the '1 + 1' photon terms because Eq. ( 140bt 
can be directly derived from Eq. ( 1321 ), assuming that the two-photon 
profile is given by independent (non-interacting) resonances, where 
the line shapes are given by the normal Voigt-profiles. 

Substituting n kpls = n p[ (v kplfi ) + A«tpi s , and using the relation 



Jcv/kTy 



[1 + /7 pl (v)]/n pl (v), Eq. flObt simplifies to 



dXj 
~d7 



2y(l+l) 



lsip 



■ A J 



kp Pel 



H X -hV fl /kTy 



+ —Xi s Aj kp p l £ a n<' 1 (Vj k p)An kvls , 
gu 



(41) 



such that upon summing over the intermediate kp resonances we 
have, 



j<kp 
Pi 



k P 



•-ftj*-. (380 



kp kp 



The partial rates in and out of the A:p-state are 



=A jkp [l+n p, (v jkp )]X J 

gkp 



(38d) 
(38e) 



At 



2y(l + l) 



,2y(l+l),< 
Vis 



AY 



+ ^X Is y A jkp p* p n pl (v;*p) An*p ls . (42) 



gu 



k=2 



The above with Eq. d34t reveals the correction term for the rate 
equations as; 



such that ■ 



S2_ = n{ e*W«V with/;/ = fiif. 

, g£x lt l L sj x ls 

When more than one two-photon channel is included, then for 
every Lyman-resonance the following needs to be computed; 



p7 = Pi 



An 



kp 



kp _ V Lkp 
1 Z_J Pi 

j 

l R l P - 



-hvtp/kTy 



3X ls R- _ V . rJ - 
n kp n kp 

where the sums run over all involved initial levels j. 



(39a) 



(39b) 



X ue - hv i l > /kT v - Xj 



+ A 



2y gj 



./Is 



AY 



gu 



j A 

A^-Z%7^" pl ^' p ) A ^ 



(43) 



Here we define A§[ = @' 2 - A^ 1+1) '*/A^ s . 

Similar to the 2s- Is two-photon channel, the correction to the 
rate equations here has two contributions. The first is related to 
AQL, which is independent of the solution to the photon distribution 
and therefore can be pre-calculated, and the second arising from the 



3.4.4 Correcting the net rates in the multi-level atom 

Equation {32) relates the population of level j with the ground state. 
The corresponding net two-photon transition rate includes the effect 



14 IChluba & Sunvaevl <2010bh proposed a varied treatment in which the 
' 1 + 1' photon terms were first taken out of the standard network of rate 
equations and then the full two-photon rate between level j and Is added, 
which at the end, is completely equivalent. 
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integral AQ\. However, in contrast to the 2s- Is two-photon chan- 
nel, in the normal rate equations, part of the latter term is already 
included. Thus the ' 1 + 1 ' photon term has to be subtracted (last 
term in brackets), where this term is calculated using the Sobolev 
approximation for An*pi s , 

A£r' in principle also arises in the normal '1 + 1' picture, when 
differentially accounting for the effect of stimulated emission in 
the CMB blackbody. However, the shape of the two-photon emis- 
sion profile is crucial, since with the normal sum of Lorentzians 
the integrand in E q. Olct would diverge for v — > Vj\ and v — » 
dChluba & Sunvad l2010b ). Furthermore, the latter two terms in 
Eq. d43t account for both, modifications in the shape of the full 
two-photon profiles, and differences in the solution of the photon 
field in comparison with the standard Sobolev approximation. 

The problem is numerical because two large terms are being 
subtracted. One way to achieve stable results is to split the range 
of integration into intervals where the mean occupation number 
in the standard Sobolev ap proximation is represented by (compare 



in the standard sobolev app 
IChluba & Sunvaevll2009bl) 



A«ipis = An 



V>') 1 



' - e - T >-*> I dv' 



An 



e 'i 



An^ [!-/>*"], (44) 



with^J p = C ip k y dv'. Outside the resonances one can simply com- 
pute each term in Eq. ( 143b separately, since there the contributions 
are small. For those intervals containing a resonance k on the other 
hand, one should compute both contributions in one integral, so that 
the main terms cancel. Clearly, the choice of the intervals is only 
motivated by the numerical precision that needs to be achieved. 
Since the Voigt-profiles have their main support inside the Doppler 
core, it is sufficient to define regions of a few Doppler width around 
the resonances. This approach suffices for our purpose. 

Alternatively, one can directly integrate the net two-photon 
production rate, Eq. l |36t , over frequency and then subtract the net 
'1 + 1" photon rate to obtain the correction. We confirmed that both 
approaches lead to the same answer. 

To capture part of the dependence of AR 2 % iU on the solu- 
tion for the populations, in numerical computations we tabulate the 
function T 2 ^ ls = AR 2 1J(X H An[) as a function of redshift, once 
we computed the solution for the photon field using the results for 
the populations of the levels obtained from a run of our effective 
multi-level recombination code. 



3.5. 1 Net rates for ns-ls and nd-ls Raman-scattering 

The net change in the number density of electrons in the level 
j 6 {ns, nd} caused by j-ls Raman scatterings is given by 



dXj 



■it**. 


f v lsc 

f 




f *< 

Jo 




f 

>yju 






-4a r 









v fls )n(v)[l + n(y - v jls )] dv 



I dv 



n(v- v flB )[i +/z(v)]dv,(45b) 

where Eq. J45bt was simply obtained from Eq. \45ai by transform- 
ing the frequency range of the second integral. 

In Eq. l !45t , ipj denotes the j-ls Raman-scattering profile, and 
the Raman-scattering coefficient is given b4^1 



f^kp j fern • 

k=n f +l %j 



(46) 



The ratio of the statistical weights is gk P /gj = 3 for the ns-states, 
and gkp/gj = 3/5 for /id-states. Equation J46l > simply reflects the 
one-photon terms and branching ratios of all the '1 + 1' photon 
routes j — > np — > Is via intermediate p-states with n > nj. 

A rigorous treatment of Eq. d46t would include the integral 
over continuum states. However, any electron reaching the con- 
tinuum would forget its history because of fast Coulomb interac- 
tions resulting in decoherence of the Raman process in the contin- 
uum. Furthermore, as mentioned above, the Lyman-continuum is 
extremely optic ally thick suc h that these channels will always can- 
cel out (see also|Hirata 2008]). Also, in numerical computations we 
only follow the evolution of the photon field up to some maximal 
frequency, v max . Therefore, in our description we are not account- 
ing for the full Raman-process connected with transitions involving 
photons with v > v max . This approximation is fully justified as the 
higher Lyman-series contribute negligible amounts to the total re- 
combination rate. Thus the sum over intermediate p-states become 
finite, without significant loss of precision. 

To simplify Eq. l |45l > we define the following quantified 



Vjis)/(v)dv 



<R\ 
K 2 



(n[l+»'])f 
<«'[!+ n] )f 



(47a) 

(47b) 
(47c) 



with n' = n(v - Vj\ s ). In the spirit of the two-photon emission and 
absorption process, we can now write 



3.5 Raman-scattering 

In our previous works IChluba & Sunvaevll2010bl. l2009j) we did 
not consider the effect of Raman-scattering on the ionization his- 
tory. However, co rrection due to this process reaches AN e /N c ~ 
0.9% at z ~ 900 dHirata|[2008t) . and hence demands careful con- 
sideration. The matrix element for this process is directly related to 
the one for the two-photon emission process by crossing-symmetry. 
In Appendix lAl we explain how to compute the Raman-scattering 
profiles, <pj(v), for the /-Is Raman process. Additional details can 
also be found in lHiratal d2008l) . where the importance of this effect 
during recombination was shown for the first time. 



1i\ 



ft, 



' n pl [1 + n pl '] \ R J- / A» n j- »P''l nR 



+ ^A«[1 +n<" ] 

R 



e -hv Jls /kT y +(An[l+n p1 ']). 
(V'fl +n pl ] 



(48a) 
(48b) 



15 We call A R ls 'coefficient' since in vacuum there is no Raman-process. 

16 Formally, the upper limit of the integral over the Raman-profiles should 
go to infinity. However, since we are following the spectrum in a finite range 
of frequencies, this introduces an upper limit, v miix < V[ sc . 
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The total j —> Is '1 + 1' Raman-scattering rate in the CMB black- 
body ambient radiation field is 



ifc=n,+l °> 



(49) 



This then leads to 

R 



= A 



Sis 



AK 



X ue - hv fls /kT ? - Xj 
{ = f <pf(v - v„ s ) [1 + «<"'] An(v) dv. 



+ — A*X ls A7?{ (50a) 

Sis 



(50b) 



Here A^ - * and A7?j are important in defining the correction to the 
rate equations (see § |3.5.4| l. Again A*'* is the thermal contribution, 
while A^j arises from non-thermal photons. 



number this equation becomes 
dAn x fl 



dz 



<T R . NnC v 2 , d> K ' r 

-£-±-^\[ Vll An{ m -ttAn<\, 
H[l + z] Vjj X 2 L J 



(53) 



where A«£ m is defined as in Eq. l |37a| l. The /-Is Raman-scattering 

]2 .R' 

cross section is given by cr , = — -£ — — , and we set d> ■' = 

4nvji<jr'*(v - V/is). Note the close similarity of this equation to 
the one-photon equation for the Lyman-series emission and ab- 
sorption channels in Eq. l !14b . Photons scattering from frequencies 
< v < vjc into the range vju < v appear as a source term. This 
is related to the fact the these photons are drawn from the CMB 
blackbody. 



3.5.2 Terms in the radiative transfer equation for 
Raman-scattering 

From Eq. ( 145 b the terms in the radiative transfer equation for the 
photon field can be obtained. However, one aspect is important to 
keep in mind: a photon that Raman-scatters off an electron in the 
/ h -state is removed from frequencies < v < v Jc . However, the 
scattered photon appears in the frequency range vj c < < V < 
Vic, and likewise for the inverse process. This description assumes 
complete redistribution of photons over the full Raman-scattering 
profiltOll during each scattering event. Therefore, the terms for the 
radiative transfer equation read 



I ^1 

c dt 



1 Wy ^ 

c dt 



Hi 

gu 



Af u N ls vJ(y)n( Vfls +v)[l + n(v)] 



(51a) 



- A%Nj tf(y) My) [1 + n(y fl , + v)] 
= A%Nj iftty - v, ls ) n(v - Vfl.) [1 + Mv)] 
^A%N ls <p*(v-v fls )Mv)[l+Mv-Vju)], (51b) 



where {pj(v) = <ffi(v)/4n. It is clear that the total integral over fre- 
quency vanishes, when adding the above two terms, showing that 
the Raman-process conserves photon number. However, the num- 
ber density of electrons in the Is and y-state is altered after each 
Raman-scattering event, according to Eq. ( 145b . 

With regards to the recombination dynamics we are not in- 
terested in the changes to the photon spectrum at low frequencies. 
Therefore, we only consider Eq. d51bb . For stimulated terms, the 
distortions can be neglected. Furthermore, one can define y£'*(v) = 
^(v) ra pl (v) [1 + n pl (v j[s + v)] ~ ^(v) n pl (v), and neglect the distor- 
tions at low frequencies, such that 



I My 

c dt 



R.Vjl S <V 



^A%Njvf(y-v Jls ) 



—A%N ls ipf(v - vj ls )fi Mv), (52) 

gls 



where // is defined by Eq. d37bb . In terms of the photon occupation 



17 We neglect corrections caused by partial redistribution in Raman- 
scattering events, but like in the case of two-photon transitions these should 
be very small. 



3.5.3 Correcting the Lyman-series emission and absorption 
terms in the radiative transfer equation 

Like in the case of two-photon emission and absorption, the reso- 
nant part of the Raman-process is already part of the ' 1 + 1 ' pho- 
ton Lyman-series transfer in Eq. < 14b . To avoid double-counting we 
simply have to correct the death probability and A/j^m of the Lyman- 
k resonance for terms that are included in the Raman-scattering 
process. For example, when using the terms for the 2s- Is Raman- 
scattering process in the radiative transfer equation, p k ^ and A«em 
for Lyman-/3, y, S, and higher will have to be corrected. 

The modified death probability can be obtained by adding 
appropriate terms to the sums of Eq. ( |39b . However, for each in- 
cluded Raman-channels one now has R 



kp 



A^[l + "P'OWL 



3.5.4 Correcting the net rates in the multi-level atom 

Like in the case of two-photon emission and absorption events, cor- 
rections to the net rates in the multi-level atom have to be defined 
to avoid double-counting. In the standard multi-level description 
of all j «-> rtp «-» Is sequences (n, < n), the contributions to the 
Raman-scattering net rate, Eq. ( 145b . takes the form 



dXj 



dXj 
~d7 



R(l + 1 



gkp 
gu 



-*ls Atpls P d n kp[s 



^lX J A kpj n p \v jk )p'Z (54a) 
gj 

gkp 



= — X u A* pi h kpls - ^Xj A kpj nP l (Vj k )pl (54b) 



R(l + 1) 



gls 



I 



dXj 
~~d7 



R(l+D 



(54c) 



As mentioned above, n max is a consequence of the finite computa- 
tional domain. The terms in Eq. ( 154b is interpreted as in the case of 
two-photon emission and absorption (see § 13. 4.4b . 

Inserting ;j* p i s = w pl (v,ipi s ) + Ah kp i s , and using the relation 



JivlkTy _ rl + „pi(y)]/„pi( v ) ; Eq. < f54bt simplifies to 

R(l+1 
ls,*p 



dXj 
~dt 



^A kpj nv\v jk )p^ 
8j 



—X u e~ hv '' lkT y 
Sis 



X, 



, g k P v At kp . - 

+ — XuA^jP^An^u, 

$ Is 



(55) 
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such that summing over the intermediate kp resonances leads to, 
gj_ 



dXj 



Is 



Vis 



+ X ls J] -A; pj ^AS ¥ , (56) 



gi 

Using Eq. J50t it is clear that the correction term for the rate equa- 
tions are, 



AR 



_/^41S 



Sj_ 

gu 



X u e 



-hvjis/kT r _ ^ 



gls 



■y gk P 



4* n v 



A;; 



(57) 



where we define AH J 2 = K{ -A^ ls + */Aj ls . Like for the two-photon 

channels (see Eq. ((43}), A??, in principle also arises in the normal 
'1 + 1' picture, where the shape of the Raman-profile ensure that 
the integrand remains finite, this time in the limit of v — » vju (see 
Appendix IA2I >. Furthermore, the latter two terms in Eq. l |57| > ac- 
count for both, modifications in the shape of the Raman profiles 
with respect to the normal sum of Lorentzians, and differences in 
the solution of the photon field with respect to the Sobolev approx- 
imation. 

Again one can compute the integrals over frequency by split- 
ting the range of integration and using Eq. d44t to model the 
"Sobolev part". In numerical calculations, we tabulate < F* <ls = 
ARj^ ls /(X u An J L ) versus redshift, to include the correction into the 
effective multi-level recombination code, and then use this to cor- 
rect the rate equations. 



4 CHANGES IN THE LYMAN-SERIES DISTORTION 
FOR DIFFERENT PHYSICAL PROCESSES 

In ^5]we discuss the changes to the free electron fraction due to the 
various physical processes under consideration. However, in order 
to understand the source of these corrections it is illustrative to first 
look at the modifications in the Lyman-series spectral distortion. 

In Fig. [TJ we present the spectral distortion at two different 
redshifts, one before the maximum of the Lyman-series emission 
(which happens at z ~ 1300 - 1400), and one just before the 
maximum of the Thomson visibility function. We include Lyman- 
resonances up to n = 8 for these computations. The solutions to the 
populations of the hydrogen levels were obtained from our imple- 
mentation of the effective 400-shell recombination code. 

The solid black line in all panels, shows our reference case, for 
which the Lyman-series is modelled using Voigt-profiles. This case 
already includes the effect of resonance scattering (for all Lyman- 
series resonances), electro n scattering, the full time-dependence 
dChluba & Sunvaevl l2009bT ) of the emission and absorption pro- 
cess, and the thermodynam ic correction factor for each resonance 
dChluba & Sunvaevl l2010bh . capturing a large part of the correc- 
tions with respect to the Sobolev treatment. In particular, the dis- 
tinction between scattering and emission/absorption events (by in- 
troducing the death probability) is important for the photon distri- 
bu tion on the blue side of the Lyman-o resonance (see discussion 
in lChluba & Sunvaevl l2009bT) . Furthermore, time-dependence and 
the thermodynamic correction factor lead to a large modification of 
the photon distribution with respect to the standard Sobolev case. 

We will now discuss the effect of the different processes on 
the shape of the Lyman-series distortion separately. 



4.1 Effect of Lyman-series scattering 

In Fig.[T] the dotted curve shows the case for which we "switched 
off" the terms for Lyman-series scattering. This line is only visible 
in the upper panels, since at high frequencies above the Lyman-o 
line it coincides with the reference case. The figure illustrates that 
partial redistribution by Lyman-series scattering is only important 
close to the Lyman-o- resonance, and on its red wing. We could, 
in principle, neglect the correction due to resonance scattering for 
Lyman-ra with n > 2, however, with our efficient PDE solver it is 
straightforward to take them into account. 

The physical reason for this behaviour is that the scattering 
probability in the Lyman-o- line is very close to unity (p^ ~ 0.999- 
0.9999), such that only in the vicinity of the Doppler core can real 
emission and absorption terms act efficiently, strongly redistribut- 
ing photons over frequency. Outside the Doppler-core, however, re- 
distribution is much slower making the effect of Doppler redistri- 
bution visible. 

For the higher Lyman-series resonance, on the other hand, the 
death probability is only about an order of magnitude smaller than 
the scattering probability, implying that far out in the wings of the 
resonance photons can be efficiently redistributed by emission and 
absorption processes. In this case, res onance scattering leads to a 
small correction (see also arguments in lAli-Hai'moud et alj20ld) . 



4.2 Two-photon emission and absorption from the excited 
states with n > 3 

Next we include the corrections due to the shapes of the rcs-ls and 
nd-ls two-photon profiles (Fig. [T] red/dashed line). By shape we 
also address modifications caused by the presence of CMB black- 
body photons. 

One can see that in comparison to the reference case this 
slightly decreases the spectral distortion between all Lyman- 
resonances, indicating that the emission/absorption opacity has 
decreased. The largest effect is seen between the Lyman-o- and 
Lyman-/? lines, as a result of the 3s-ls and 3d-ls two-photon emis- 
sion and absorpti on process. This result is in agreement with our 
earlier treatment dChluba & Sunvaevll2010d) . where it was demon- 
strated that the shape of the 3s- Is and 3d- Is two-photon profiles 
leads to a slight acceleration of recombination, which however, is 
less important than the corrections arising from the thermodynamic 
correction factor and time-dependence, individually. 

We tried to identify the main source of the modifications above 
the Lyman-/? resonance in more detail. In the full two-photon pic- 
ture, the 3s- Is and 3d- Is two-photon emission and absorption chan- 
nels only act on photons with v < v 31 . However, when neglecting 
the modifications to the shapes of the two-photon profiles, a large 
part of the opacity above the Lyman-/? line (incorrectly) comes 
from the 3s- Is and 3d- Is '1 + 1' photon channel, which involves 
the Lyman-o resonance and is modelled by a normal Voigt profile. 
It turns out that only for the 3s- Is and 3d- Is two-photon process 
does the exact shape of the two-photon profile really matter. Above 
the Lyman-/? line the small correction with respect to the solid line 
is practically captured by truncating the Voigt profiles (in particu- 
lar the one for Lyman-o), such that the energy is conserved (e.g. 
photons emitted or absorbed in a 3s- Is and 3d- Is '1 + 1' photon 
process can only have energies v < v 3] , and so on). This illustrates 
how important the shape of the line profiles is when going far into 
the damping wings of the resonances. 
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Figure 1. Solution for the Lyman-series distortion at z = 1490 (left panels) and z = 1190 (right panels) for different combinations of physical processes 
(for details see !j4}. In all cases we include the effect of electron scattering. We also marked the positions of the Lyman-series resonances with the vertical 
dashed-dot-dotted lines. The effect of partial frequency redistribution is only important close to the Lyman-a line center, so that the dotted line is only visible 
in the upper panels. A movie on the time-evolution of the Lyman-series distortion can be found at www.Chluba.de/Lyman-series-movie. 



4.3 Importance of Raman scattering 

In addition to the two-photon corrections, we ran cases that also in- 
clude the full Raman scattering treatment (Fig.Q] blue/dash-dotted 
line). One can see that the Raman process led to an enhancement 
of the spectral distortion between the Lyman-o- and Lyman-/? reso- 
nance, while in all the other cases the spectral distortion decreases 
between the resonances. Thus one expects an increased blue-wing 
feedback correction and hence a delay of hydrogen recombina- 
tion from Lyman-a. On the other hand, these additional red-wing 
Lyman-/? photons were created in a 2s- Is Raman event, such that at 
earlier times an acceleration is expected. This si mple picture i s in 
agreement with earlier discussions of this process jHiratall2008h . 

In the case considered, the main source of the difference in the 
Lyman-series distortion comes from the 2s- Is Raman treatment. 
Neglecting the Raman-corrections to the higher ns-ls and nd-ls 
channels does no affect the shape of the distortion noticeably. This 
is one of the reasons why the Raman process need to be included 
only for the first few levels. 

Since in the case of Raman scattering, the 2s — > Is scattering 
profile is given by $'*(y) ~ ^(v)n pl (v) (see Sect. I3.5l l. one ex- 
pects two sources of corrections: (i) due to the difference of tpj(v) 
with respect to a sum of Voigt profiles with appropriate weights, 
and (ii) the factor n v \v). In the normal '1 + 1' photon picture this 
factor would not appear differentially, but instead directly for each 



resonance frequency. It turns out that both part of the correction are 
important for the 2s- Is Raman treatment. 

We note, that the spectral distortion at z = 1 190 in the fu ll 
treatment lo oks very simi lar to the curve given in iHiratal d2008h . 
However, in lHiratal d2008l) also the CMB blackbody spectrum was 
added, and n v instead of x 3 An v was plotted, which makes a direct 
comparison more difficult. 



5 CHANGES TO THE FREE ELECTRON FRACTION 
FOR DIFFERENT PHYSICAL PROCESSES 

In this section we present our analysis of the different corrections to 
the standard recombination calculation. We focus on hydrogen, and 
model t he helium recombination d ynamics using the description 
given in IChluba & Sunvaevl fcoioj) . including the first five shells 
with full feedback. With the current version of our effective multi- 
level recombination code we are able to account for all important 
corrections to the recombination dynamics of hydrogen. We show 
a direct comparison with previous results and find very good agree- 
ment. All figures in which we compared the output of our recom- 
bination code with Recfast we used Recfast vl.4.2 JWong & Scottl 
|2007|) . but excluded the corrections to the helium recombination 
history in Recfast and removed the switches in the Recfast ODE 
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Figure 2. Dependence of the modifications to the recombination dynamics 
on the number of included shells. The results of our effective multi-level 
recombination code were directly compared with the output from Recfast. 

system (see Fendt et al. 2009, for details). The cumulative effect on 
the ionization history is presented in Fig. [6] 

5.1 Results from our effective multi-level code 

In Fig.|2]we show the changes in the recombination dynamics with 
the number of shells that were included into the computation of 
the effective rates. This figure confirms that our implementation 
of the effective multi-level approach yiel ds corrections that a re in 
agreement with our earlier computations IChluba et alj l l2010h . We 
find that the correction converges down to z ~ 200 when includin g 
~ 300 - 400 shells, as already explained in lChluba et"aTI J2010l) . 
We also directly compared with our full multi-level recombination 
code and found the difference to be smaller than AN C /N S ~ 10" 5 . 

Collisional processes are still able to change th e low redshift 
behav iour at the ~ 0.1% level in this redshift range dChluba et al.1 
however, we defer a detailed analysis on the importance of 
this effect to a future work. 



5.2 The reference case 

In Fig. [3] we present a compilation of different corrections to the 
ionization history that are included into our reference case. For 
this we internally compared the outputs of our recombination code 
when switching on and off different processes. We computed the 
solution to the photon transfer problem including the Lyman-series 
up to n = 8, with Lyman-# (n = 9) on the upper boundary of 
the frequency grid. In the Lyman-series transfer we did not in- 
clude the corrections to the profiles of the emission and absorption 
processes arising from two-photon and Raman events, i.e. we de- 
scribed Lyman-fc emission and absorption using Eq. dl4t . However, 
in our full reference case, Lyman-fc resonance and electron scat- 
tering, as well as 2s- Is two-photon emission and absorption were 
included (see Sect.|4]for additional comments). 

To account for all the corrections to the rate equations in the 
effective multi-level recombination code, we ran the obtained solu- 
tion for the photon distribution through the modules that also allow 
us to take the two-photon and Raman scattering corrections into ac- 
count (see explanations in Sect. [3}- However, we replaced the full 
profiles of the channels with the normal Voigt-profiles. 



Figure 3. Corrections that are included into the reference case. 

The cumulative difference with respect to the output of our ef- 
fective multi-level recombination code which does not include any 
of the radiative transfer corrections is shown in Fig. [3] In total we 
find a delay in recombination by AN e /N e ~ 0.4% at z ~ 930, and an 

early acceleration by AN C /N C 1.0% at z ~ 1270. The reference 

case therefore captures a significant part of the total correction with 
respect to Recfast (see Sect. |5.5| for details). 

5.2.1 The 2s-ls two-photon correction 

Figure [3] shows the total correction due to changes in the 2s- Is 
two-photon channel. We only modified the 2s- Is two-photon and 
Ly-or net rate in our effective multi-level recombination code using 
Eq. d26l >. but did not alter any of the other rates. Also we switched 
off line-diffusion. 

We find a delay of recombination by AN S /N C ~ 0.83% at z ~ 
990, which is slightly (by AN e /N ? ~ 0.18%) larger than in earlier 
computations of this process (e.g. lFendt et alj|2009l) . 

There are two main reasons for this difference; (i) because in 
the reference case we include the emission and absorption in the 
2s- Is two-photon channel, the self-feedback of photons emitted by 
2s- Is transitions on the ls-2s two-photon channel is accounted for, 
which leads to an additional delay of AN C /N e ~ 0.08% and (ii) 
the remaining deceleration by AN S /N C ~ 0.1% is just caused by 
normal absorption of 2s- Is photons by Lyman-a (without the aid 
of line-diffusion). 

5.2.2 Feedback to Lyman-a and the diffusion correction 

We now consider ed the feedback correction to Ly-cr. Like in our 
earlier treatment l lChluba & Sunvaevl |2007t) we find ANJN, ~ 
0.21% at z ~ 1 100. We computed this correction from our radiative 
transfer code by modifying the Lyman-a escape probability, with 
resonance scattering switched off. We also left the rate equations 
for the higher Lyman-series resonances unaffected, in order to not 
reflect the full Lyman-series fee dback correction, which a mounts 
to ANJN, ~ 0.26% at z ~ 1 100 l lChluba & Sunvaevll2010al) . 

In Figure [3] the correction due to Lyman-a diffusion alone is 
depicted. Again this was computed as a correction to the Lyman- 
a resonance only. We find an acceleration by AN e /N t ~ -0.44% 
at z ~ 900. This is slightly smaller than in our earlier computation 
dChluba & Sunvaevl2009d) . The reason is simply that there we only 
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Figure 4. Two-photon corrections from highly excited levels (n > 2) and 
their convergence with n. 



included 3 shells into our computation. However, when including 
more than 5-10 shells the diffusion correction becomes slightly 
smaller, reducing to the curve presented in Fig. [3] To check the 
precision of our PDE-solver, we recomputed the curve for the 3- 
shells case and confirmed our earlier result. 



5.2.3 The correction due to time-dependence and 
thermodynamic factor 

The final correction that is taken into account by the computa- 
tions in the reference case is caused by the time-dependence of the 
Lyman-o- emission process and the thermodynamic correction fac- 
tor (se e Fig. [3). The origin of thes e ter ms was first explained in de - 
tail bv lChluba & Sunvaevl J2009bh and lChluba & Sunvaevl feOlObl) . 
The net effect is an acceleration in recombination by AN e /N s ~ 
-1.28% at z ~ 1200. This result is in excellent agreement with 
the curves presented in IChluba & Sunvaevl fcfJlObh. Fig. 18 therein. 
Note that in the Figure o f lChluba & Sunvaevl d20 lObh also the 3s- Is 
and 3d- Is two-photon profile correction was included. 



5.2.4 Correction from the higher Lyman-n 

For the solid black line in Fig.[3]all Lyman-series corrections were 
included. However, so far we have just discussed the corrections to 
the 2s- Is two-photon and Lyman-Q' channel. We found a cumula- 
tive acceleration by A/V c //V c 0.06% at z ~ 1210 as result of the 

higher Lyman-series. This correction includes all feedback correc- 
tions among the higher levels, time-dependence, and the thermody- 
namic factors. Since the additional modification is small, this shows 
that just a detailed treatment of 2s- Is two-photon and Lyman-a cor- 
rections already gives a very good approximation to the total cor- 
rection in the reference case. 



Figure 5. Raman scattering corrections from excited levels 

field. For two-photon processes from ns and nd-states with n > 2 
the former dominates. 

In Fig. [4] we present the changes to the free electron fraction 
when the two-photon corrections up to the 8s- Is and 8d-ls two- 
photon channel are included. Two-photon processes lead to a total 

acceleration of recombination by &N e /N e 0.46% at z ~ 1120. 

The main contribution comes from the 3s- Is and 3d- Is two-photon 
process, while the higher levels only add A/V c /7V e ~ -0.08% at 
z ~ 1200. The correction practically converge when accounting for 
the two-photon terms up to 5s-ls and 5d-ls. Also the result for the 
3s-ls and 3d-ls two-pho ton process compares extre mely well with 
our earlier computation dChluba & Sunvaevll2010br) . We conclude 
that for practical purposes it is sufficient to include the two-photon 
corrections for all ns and nd states up to n ~ 4 - 5. 

5.4 Corrections caused by Raman processes 

The final process we discuss is the effect of Raman scattering, 
which was investigated also by iHiratal fa)08l) . The result of our 
computation is shown in Fig. [5] confirming that Raman scattering 
leads to a delay of recombination by AN C /N e ~ 0.9% at z ~ 920. 
This r esult is in very good agreement with the analysis of IHiratal 
d2008l) . We found that the correction is dominated by the 2s- Is Ra- 
man process. Higher level Raman scattering lead to a small ad- 
ditional modification, which for practical purposes could be ne- 
glected. We recommend including the Raman-corrections for the 
first three shells. 

As mentioned above, the Raman correction has two separate 
contribution: one from the feedback of Lyman-/? photons on the 
Lyman-Q' resonance, which leads to a delay of recombination a 
low redshifts, and a second arising because of the accelerated 2s- 
ls scattering. The delay and acceleration need to be out-of-phase 
in redshift to create a net effect. We find that their individual con- 
tributions have amplitudes comparable to ~ 2%, however, partial 
cancelation makes them smaller. 



5.3 Two-photon corrections from levels with n > 2 

In this Section we discuss the correction caused by the modifica- 
tions in the emission and absorption profiles of the ns-ls and /id- 
Is two-photon channels. These corrections are due : (i) quantum- 
mechanical modifications to the shapes of the line-profiles, and (ii) 
stimulated two-photon emission in the CMB blackbody radiation 



5.5 Total correction with respect to Recfast 

In Fig. [6] we show the cumulative correction to the ionization his- 
tory caused by all the processes included into our present recombi- 
nation code. The changes during hydrogen recombination found 
here are very similar to those presented in iRubino-Martfn et "all 
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Figure 6. Total correction to the ionization history. We compared the out- 
put of our effective multi-level recombination code with Recfast vl.4.2 
IWong & Scottll2007t) . However, we switched all helium flags in Recfast 
off. 



and to optimize the implementation of the recombination code, 
so that runtimes of seconds can be accomplished, incorporat- 
ing all the important physical processes without requiring any 
fudge factors. Our final version of CosmoRec will be available at 
www.Chluba.de/CosmoRec. 
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The only major difference is visible at low redshifts, since 
their analysis was based on the results from a 100-shell hydrogen 
recombination model and hence the low reds hift freeze-ou t tail w as 
overestimated (see Fig. [3] a nd comments in Chluba et al. 1 2010h). 



However, as argued earlie r dFendt et alfcOOsj ; iRubino-Martm et al .1 
l201fj ; IChluba et alJuOlOf) . we expect that this additional modifica- 
tion does not affect the conclusions of their work at a significant 
level. In particular, they showed that the cumulative effect of all 
published recombination corrections could lead to biases in the val- 
ues of fib h 2 and n$, that reach l.lcr and 2.3cr, respectively. 

For the analysis of CMB data from the Planck Surveyor these cor- 
rections have to be taken into account carefully when answering 
queries about different models of inflation. 



6 CONCLUSIONS 

In this paper we complete our analysis on the importance of two- 
photon transitions and Raman scattering during the cosmologi- 
cal recombination epoch. We explicitly solve the radiative trans- 
fer equation for the Hi Lyman-series transport, including all im- 
portant processes (e.g. resonance scattering, full time-dependence), 
extending our former treatment to account for Raman scattering, as 
well as two-photon transition from highly excited levels with n > 3. 
Our computations are performed using an effective multi-level ap- 
proach for hydrogen to accelerate the recombination calculation, 
that without optimization achieves runtimes of ~ 1 — 2 minutes. 

We find that 2s- Is Raman scattering leads to a small initial 
acceleration of recombination at high redshifts, which then turns 
into a deceleration of AA C /A e ~ 0.9% at z ~ 920. ras-ls and 
rad-ls Raman processes from levels with n > 2 only result in a 
small additional correction. Two-photon transitions from ns and nd- 
states with n > 3 accelerate hydrogen recombination by additional 

AAe/A'e 0.08% at z ~ 1200. For practical purposes one only has 

to include the two-photon corrections for the first ~ 4 - 5 shells. 

This work carves a path towards a new cosmological re- 
combination code, CosmoRec, that supersedes the physical model 
included in Recfast and can be used in the analysis of future 
CMB data, e.g. from the Planck Surveyor, Act, Spt, and CmbPol. 
The final step will be to perform a detailed code comparison, 



APPENDIX A: COMPUTING THE TWO-PHOTON 
EMISSION AND RAMAN PROFILES 

Al Two-photon emission profiles 

We compute the two-ph oton decay profiles a c cording to 
IChluba & Sunvaevl J2009bl) and IChluba & Sunvaevl feOlObh . In 
their treatment the infinite sum over intermediate /7p-states is split 
up into levels with principal quantum numbers n > ;z, and n < n h 
where n, is the principal quantum number of the initial state. This 
makes the sum over resonances (in the case of 3s and 3d only one) 
finite and allows to use interpolation or fitting formulae for the re- 
maining contributions to the total matrix element from the infinite 
sum. 

We tabulate the non-resonant parts of the two-photon ma- 
trix elements prior to the computation. The resonances are an- 
alytical ly added afterwards. As explained in IChluba & Sunvaevl 
d2010bh . close to the resonances, the motion of the atoms becomes 
important, leading to a broadening of the two-photon profiles. To 
include the effect of motion on the shape of the lines close to the 
Doppler core, we take the ratio, p 2y = ip 2y /<pf A , of the vacuum 
two-photon profiles to the sum of Lorentzians, and tabulate it on 
the computational grid in frequency. For every evaluation in time, 
we first compute the Voigt-profiles for the resonances of interest 
and then sum these with their respective weights to obtain the total 
Voigt-profile of the resonance. These are then multiplied by pf' to 
obtain an approximation for the two-photon profile in the lab frame. 

This procedure also allow us to include the changes in the total 
width, r„ p , of the resonances with redshift. In vacuum this width is 
related to the total decay rate, however, with the CMB this rate can 
change at the level of ~ 10% for the Lyman-n line when n > 2. 



A2 Raman profiles 

A ns-ls and nd-ls Raman process has form H* + y — > H + y', 
where H* denotes a neutral hydrogen atom in an excited ns/nd- 
state. In contrast to the two-photon emission process, H* — > H + 
y + Y, the Raman process only works when photons are available. 
In vacuum no Raman scattering events occur. The ns-ls and nd-ls 
Raman scattering matrix elements are related to the ns-ls and rcd-ls 
two-photon matrix elements by crossing-symmetry. The energies of 



© 0000 RAS, MNRAS 000, 000-000 



16 Chluba and Thomas 




V / V 2 



Figure Al. Frequency dependence of the 2s- Is Raman profile. This profile 
has to be interpreted as a scattering cross section, where the electron is in 
the excited state. The positions of the Balmer-cc, fi, 7, and S resonances are 
marked. For comparison the cross section computed as a sum of Lorentzians 
are shown. Furthermore, we also show the ratio of these two profiles, indi- 
cating that close to the resonances the profile becomes Lorentzian, with 
small corrections. 



the incoming photon, 7, and the outgoing photon, 7', are given by 
v + v jH = v ', where hvj is is the excitation energy of the initial state 
with respect to the ground state. 

Ra man scattering profiles can be derived using the formulae 
given m IChluba & SunvaevlfcOlObl) . The main difference are: (i) in 
the functions f„(y) and h„(y) (see Eq. (8c) and (lOd) in their paper) 
y = v/Vji s has to be replaced with —y; (ii) the pre-factor y 3 (l - 
y) 3 needs to be substituted by y 3 (l + y) 3 ; and (iii) the resonances 
now appear for intermediate zip-states with 11 > where iij is the 
principle quantum number of the initial state. 

For a given computational frequency grid only a finite number 
of resonances appear, say for rij < n < n rcs . Like in computations 
of the two-photon emission profiles, one can therefore split the infi- 
nite sum over intermediate p-states (including the continuum), into 
resonant and non-resonant contributions. The non-resonant contri- 
butions come from n < nj and n rcs < n, where the non-resonant 
matrix element scales like 1/y for y — > and (l/[n 2 . - 1] - y)" 1 
towards the ionization threshold, v — > Vj c . One can therefore tab- 
ulate M m y [l/[rrj — 1] - y] on a grid and add the finite number of 
remaining resonances analytically. The pole displacements arising 
from the finite lifet i me of t he intermediate p-state, as mentioned in 
IChluba & Sunvaevl ( feOlObl) . has to be included. 

In Fig. lAll and lA2l we present some examples of Raman scat- 
tering profiles. The electron is assumed to be in the excited state, 
so that a low frequency photon can Raman scatter off the atom. In 
the recombination problem, these photons will be drawn from the 
CMB blackbody, as spectral distortions below the Balmer contin- 
uum can be neglected. Figures [AT] and fA2] also show the ratios of 
the Raman profiles with respect to the sum of Lorentzians. Close 
to the resonances all these ratios are extremely close to unity, as 
expected. We use this ratio to include the effect of motion of the 
atom on the shape of the resonances close to the Doppler core, as 
explained in the previous section. 
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Figure A2. Same as Fig. lAll but for the 3s-ls and 3d-ls Raman profiles. 
We marked the positions of the Paschen-ff, /}, and 7 resonances. 



APPENDIX B: PDE-SOLVER 

For this work we implement our own partial differential equation 
(PDE) solver in order to fine tune performance and precision. Com- 
paring with previous results obtained using the NAG library con- 
firms the precision of our own implementation. 

The PDE describing the radiative transfer problem during re- 
combination is of the parabolic type. It is desirable to use an im- 
plicit or semi-implicit numerical scheme, to avoid strong limitation 
on the step size imposed by stability. Several numerical algorithms 
for this type of problems have been discussed, e.g. Crank- Nicolson 
method (see lAntij2002h . 

For the recombination problem it is beneficial to use a non- 
uniform grid in frequency as in the vicinity of the resonances one 
needs a resolution of A v/v ~ 10~ 7 - 10~ 5 , while a much coarser grid 
can be introduced outside this zone. However, this implies that the 
spatial discretization that is normally used in the Crank-Nicolson 
method is only accurate to first order in the grid spacing. We there- 
fore decided to implement a second order scheme in which the first 
and second derivatives of the occupation number n x with respect to 
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frequency can be written as 

^ =Y J «i(x)n Xi (Bl) 

f 

d 2 n x v-i 

-^=Y,Hx)n Xi (B2) 

where the sums run over five grid-points in the neighbourhood 
of x. The coefficients /c,(x) and A,-(jc) can be easily derived using 
Lagra nge interpolation formulae (e.g. see lAbramowitz & Stegunl 
1 1972b . These coefficients can then be precomputed and stored once 
the grid is chosen. The PDE appearing in the diffusion problem can 
thus be written as matrix equation, 

Hi.i», /'■ (B3) 

where the matrix, By, is bandec0 with four off-diagonal elements. 
Such system can be easily solved with 0{M) operations, where M 
denotes the number of grid-points. 

For each resonance we typically needed ~ 10 3 points in fre- 
quency. Increasing this number to ~ 10 4 per resonance did not 
make a notable difference for the final correction to the ioniza- 
tion history. Our typical step size in redshift was Az ~ 1, but we 
also tried a ten times smaller time-step, without finding any signif- 
icant modification in the solution. Our tests also showed that even 
Az ~ 10 should be sufficient for detailed computations of the ion- 
ization history. 

Furthermore, we tried fully implicit and semi-implicit 
schemes (^-method with G > 0.5), finding good performance for 
9 ~ 0.6. We also experimented with the distribution of grid points, 
and found that it is important to sample the Doppler core well. 
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